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Abstract 

We study the phase mixing and dissipation of a packet of standing shear Alfven waves lo- 
calized in a region with non-uniform Alfven background velocity. We investigate the validity 
of the exponential damping law in time, exp(— At 3 ), presented by Heyvaerts & Priest (1983) 
for different ranges of Lundquist, S, and Reynolds, R, numbers. Our numerical results shows 
that it is valid for (R, S) > 10 7 . 
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1 Introduction 



Since the discovery of the hot solar corona by Edlen (1943), different theories of coronal heating 
have been put forward and debated. Heyvaerts & Priest (1983), hereafter HP83, were first to 
suggest that the phase-mixing of Alfven waves in coronal plasmas could be a primary mechanism 
in coronal heating. They showed that the phase-mixing occurs due to inhomogeneity of the 
local Alfven phase speed across the background magnetic field. HP83 analytically showed that 
in both the strong phase-mixing limit and the weak damping approximation, the amplitude of 
standing Alfven waves decays with time as exp(— At 3 ) where A is a function of the coordinate 
corresponding to the inhomogeneous direction (x, in this paper). Since then, much analytical and 
numerical work has been done on the subject. Nocera et al. (1984) studied the phase mixing 
of propagating Alfven waves in an inhomogeneous medium. They pointed out if transverse 
gradients are smeared out as soon as they are formed, this yields to weak phase mixing where 
damping laws differ from solutions of HP83. Parker (1991) investigated the effect of a density 
and/or temperature gradient in the direction of vibration of a transverse Alfven wave. The result 
was a strong coupling of the waves on different lines of force, producing a coordinated mode 
that was not subject to simple phase mixing. Hood et al.( 1997a) derived a self similar solution 
of the Alfven wave phase- mixing equations for heating of coronal holes. They showed that the 
damping of the waves with height follows the scaling predicted by HP83 at low heights, before 
switching to an algebraic decay at large heights. Hood et al. (1997b) obtained a simple, self 
similar solution for the heating of coronal loops by phase mixing. They showed the HP83 model 
still does work well in a certain class of coronal loops and the phase mixing can supply heating 
at large Lundquist number at timescales shorter than or comparable with the radiative cooling 
timescale. Nakariakov et al. (1997) considered the nonlinear excitation of fast magnetosonic 
waves by phase mixing Alfven waves in a cold plasma with a smooth inhomogeneity of density 
across a uniform magnetic field. They suggested this nonlinear process as a possible mechanism 
of indirect plasma heating by phase mixing through the excitation of fast waves. But Botha et 
al. (2000) showed that the nonlinear generation of fast modes by Alfven waves has little effect 
on classical phase mixing. De Moortel et al. (1999) elaborated the effect of density stratification 
on phase-mixing. They remarked that when the inhomogeneity in the horizontal direction in 
the plasma is sufficiently large, so the phase mixing is strong, stratification is unimportant. In 
the other words due to the rapid phase mixing, energy can be dissipated before the effects of 
stratification build up. De Moortel et al. (2000) studied the combined effect of a gravitationally 
stratified density and a radially diverging background magnetic field on phase mixing of Alfven 
waves. They found that: i) The efficiency of phase mixing depends strongly on the particular 
geometry of the configuration, ii) Depending on the value of the scale height the wave amplitudes 
can damp either slower or faster than in the uniform non-diverging model. 

Hood et al. (2002) showed that the amplitude of single pulse and bipolar pulse traveling in 
the z direction, contrary to infinite wavetrain, have slower algebraic damping of the form t~ 3 / 2 
and t -3 , respectively, rather than exponential in time. Tsiklauri et al. (2003) cleared that the 
decay rate of the Alfvenic part of a compressible 3D MHD pulse is affected linearly by the degree 
of localization of the pulse in the homogeneous transverse direction, but the dynamic of Alfven 
waves can still be obtained from the previous 2.5D models, e.g., Hood et al. (2002). Smith et al. 
(2007) found that in presence of the both density stratification and magnetic field divergence, 
the enhanced phase mixing mechanism can dissipate Alfven waves at heights less than half that 
was predicted by the previous analytical solutions. They stated that if phase- mixing takes place 
in strongly divergent magnetic fields, it is not necessary to invoke anomalous viscosity in corona. 

In the present work, we study the phase mixing of a packet of standing shear Alfven waves 
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in presence of the both viscous and resistive dissipations. To do this, we numerically solve the 
linearized MHD equations and obtain the damping time of the oscillations. Our aim is to test 
the validity of HP83's damping law for different ranges of the Reynolds and Lundquist numbers. 
This paper is organized as follows. In Sect. 2, we introduce the basic equations of motion and 
introduce the model. In Sect. 3, the numerical results are reported, while the conclusions are 
given in Sect. 4. 

2 Equations of motion 

The linearized MHD equations for a zero-/3 plasma are: 

^ = -^{(V x <5B) x B + (V x B ) x SB} + ^-V 2 5v, (1) 
at 47rpo Po 



^ = V x (<5v x Bo) + ^V 2 5B, (2) 

where (5v and (5B are the Eulerian perturbations in the velocity and magnetic fields; po, a, 
7] and c are the mass density, the electrical conductivity, the viscosity and the speed of light, 
respectively. 

The simplifying assumptions are: 

• under coronal conditions gas pressure is negligible (zero-/3); 

• the equilibrium density profile is po = po(x); 

• there is a constant magnetic field along the z axis, Bo = Bqz; 

• there is no initial steady flow inside or outside of the tube; 

• the viscous and resistive coefficients, r\ and a respectively, are constants. 

To solve Eqs. (l)-(2), following HP83, we neglect the variations in y direction, = 0, and 
will further assume that the velocity perturbs in y direction. So we choose a solution for the 
velocity perturbation as 

5v(x, z, t) = 5v y (x, t) sin(kz)y, k = mr/L, (3) 

where L is the length of the loop and n = (1, 2, 3, ...) is the wave number in z direction, respec- 
tively. Here the waves are standing because of boundary conditions 5v(x, 0, t) = 5v(x, L, i) = 0. 
Note that HP83 also supposed that the loop is bounded from above an below by boundaries 
at altitude z = and z = L. It is convenient to work with dimensionless variables x = x/a, 
z = z/a, i = t/rA, Po(x) = Po(x)/poo, 5v = 5v/va and <5B = 5~B/Bq. Where a is a typical 
length scale of density inhomogeneity across the field (i.e. loop radius) and ta = cl/va is a 
time scale for an Alfven wave to propagate along the inhomogeneity direction, poo and va 
are the plasma density and Alfven speed at x = 0, respectively. Finally, Eqs. (1) and (2) in 
dimensionless form, dropping the 'bars' for convenience, become 

85v v 2 , .dSB v 1 2r , A . 

-dt = vAx) ^r + r v **> (4) 
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where va(x) = , = = is the dimensionless form of the Alfven speed. Also the Reynolds 

VP0W//900 

number, 

is the ratio of the viscous time-scale to the Alfven crossing time, and the Lundquist number, 

\ & > I \v Ao J 

is the ratio of the resistive time-scale to the Alfven crossing time. Removing 5B y from Eqs. (4) 
to (5), and keeping only the first order-terms in 1/R and l/S gives 



d 2 5v 



y- + k\\{x)5v y = (I + V 2 ^ + r(x,t), (6) 



at 2 AW y Vi? sJ at 

where 

T(x,t ) = - 6(^ -2U -4(— 7T -« T' 7 

where prime indicates a derivative with respect to x. It is obvious that the term T(x,t) in Eq. 
(6) becomes important for high magnetic diffusion plasmas (low S) and in the regions where the 
Alfven speed has a large gradient. For more accuracy, in contrast with HP83, we keep this term 
in our numerical simulations. 

Since the dissipation rate is a function of x, for calculating the overall damping time, it is 
suitable to calculate the dimensionless total energy (kinetic energy plus magnetic energy) of the 
packet per unit of length in the y direction as 

E tot {t)= [p {x)5v y '(x,t) + 5B y l (x,t)\dx, (8) 

J 

where E tot (t) = ^ L E tot (t) and 5B y (x,t) is calculated from Eq. (5). 

We suppose a functional form of dimensionless Alfven speed and a Gaussian form of a 
localized packet of standing Alfven waves around x = 1 as 

v A (x) = 2 + tanh[a(x - 1)], (9) 



5v(x, z, t = 0) = exp 



(10) 



where parameter a controls the size of inhomogeneity and d is width of the packet. For a = 2 
and d = 0.1, the Alfven speed profile and shape of the initial wave packet given by Eqs. (9) and 
(10) are plotted in Fig. 1, respectively. 

Substituting Eq. (9) in uj(x) = kvA(x) gives the dimensionless average period of oscillation 

as 



avg 



P(x)dx 



TT 



dx 



k Jo 2 + tanh[a(x - 1)] ' 



(11) 
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3 Numerical Results 



As typical parameters for a coronal loop, we assume L = 10 5 km, a = 10 3 km, Bq = 100 G, 
and poo = 2 x 10 -14 gr cm -3 . For such a loop, one finds va = 2000 km s -1 . Here the loop 
parameters coincide with the TRACE observations (see Aschwanden et al. 2002; Verwichte et al. 
2004). We use a finite difference method to solve Eq. (6), numerically. The evolution of a packet 
of fundamental standing Alfven modes is calculated in the range of < x < 2. To include the 
dynamical effect of the exterior region, we let the wave packet to evolve up to x = 2. We suppose 
that the wave packet never reach at the x boundaries. Hence to avoid any contamination of the 
solution by the change of boundary values, we fix the boundary conditions. This restricts the 
time of simulation, but it is still possible to reach the strong phase mixing limit. We choose the 
boundary and initial conditions as 

5v y (x = 0, t) = 5v y (x = 2, t) = 0, (12) 
5v y (x,t = 0) = e - 5 °( a; - 1 ) 2 , (13) 
d5v y (x, t) 



= 0. (14) 

t=o 



dt 

There is an upper limit for the time of simulation because we can simulate the evolution 
until any excitement near the x boundaries could be occurred. The truncation error of numerical 
results is A = 0(At 3 ) + 0(Ax 4 ). We should be aware of choosing suitable spatial step size Ax, 
because in the limit of strong phase mixing, large gradients in the x direction are made, so the 
smaller Ax is needed. 

Fig. 2 shows contour plots of Sv y (x, t) in the x — t plane for two different cases with R = S = 
10 4 (a) and R = S = 10 8 (b). The white and black colors represent positive and negative values 
of Sv y (x,t), respectively. Fig. 2 clears that the defocusing of the packet in the case (a) is large 
but not in the case (b). This is because of coupling of oscillations in neighboring field lines due 
to presence of damping terms in the right hand side of Eq. (6). Fig. 3 presents the cross-section 
cuts along x = 1, x = 0.6 and x = 1.4 for the case (a) with a = 2 and d = 0.1. It illustrates 
that as central regions of the packet decay with time, the neighboring oscillations in the regions 
with smaller amplitudes, are excited and finally are damped by phase mixing. This means that 
the packet defocuses along the x direction which is illustrated in Fig. 4. Fig. 5 shows the time 
evolution of the kinetic energy, magnetic energy and total energy of the packet. Fig. 5 reveals 
that both the kinetic and magnetic energies of the packet oscillate with time sharply at initial 
stage of the evolution and then smoothly damped. 

Figs. 6-8 show evolution of the packet for the case (b) with a = 2 and d = 0.1. Figs. 6-8 
in comparing with Figs. 3-5 present that in high Reynolds and Lundquist numbers, i.e. weak 
damping, the wave packet is damped in developed stage of phase mixing and it's defocousing is 
negligible. 

From Eq. (8) for R = S = 10 4 , a = 2 and d = 0.1, we obtain r^am = 79.1 s. From Eq. (11) 
for q = 2, the average period of the fundamental mode, k = ^j, is obtained as P^ 2 = 57.9 s. 
Therefore the ratio of the damping time to the average period, Td am /P^ 2 f° r the fundamental 
mode is 1.4. From Eq. (8) for R = S = 10 8 , a = 2 and d = 0.1, r da m = 1702.6 s and 
Tdam/-P a vg 2 - 29.4. For R = S = 10 s , if we set a = 4 and d = 0.05 then P a °= 4 = 62.1 s, 
Tdam = 1072.7 s and Td am / P^Jg 4 — 17.3. These strong damping times are in agreement with 
the results observed by Nakariakov et al. (1999) and Wang & Solanki (2004) deduced from the 
observation of TRACE. 
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To test the validity of HP83's damping law in the both weak damping and strong phase 
mixing limit, we fit the functional form exp(— At B ) on the envelope of Sv y (x,t) at x = 1. Note 
that B = 3 and A{x) = ^vk 2 v£(x) in HP83, where v = av a (R _1 + S^ 1 ). The numerical results 
obtained for A and B for different ranges of the Reynolds and Lundquist numbers are tabulated 
in Table 1. It shows that for R = S = 10 4 — 10 10 , the numerical values of B converge to its 
analytical value but there is one to four order of magnitude difference between the numerical 
and analytical values of A. This returns to keeping the term T(x, t) in Eq. (6) which has been 
missed in HP83. Table 1 also shows that for S = 10 12 and R = 10 7 — 10 9 , the contribution of 
T(x, t) becomes negligible in Eq. (6) and the numerical values of A and B converge to their 
corresponding analytical values in HP83's damping law. Finally one can conclude that the 
exponential damping law in time of HP83, exp(— At 3 ), is valid for (R, S) > 10 7 . 

4 Conclusions 

Phase mixing of a packet of standing Alfvenic pulses in fundamental mode is studied. Using a 
finite difference method, the linearized MHD equations for a zero-/3 plasma are solved, numeri- 
cally. The damping times of oscillations in presence of the both viscous and resistive dissipations 
are calculated, numerically. They are in good agreement with the TRACE observations. The 
exponential damping law in time of HP83, exp(— At 3 ), for the different ranges of the Reynolds 
and the Lundquist numbers are examined. Our numerical results shows that it is valid for 
(R,S) > 10 7 . 
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(a) (b) 




X 



Figure 1: (a) The profile of background Alfven speed (solid curve) and the initial amplitude of 
velocity perturbations (dash-dotted curve) at z = 50 as functions of x. (b) 3D view of the packet 
of standing Alfvenic pulses in fundamental mode (n = 1). Auxiliary parameters are: a = 2, 
d = 0.1, and L = 100a. 
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(a) (b) (c) 
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Figure 3: Cross section cuts of Sv y along (a) center of the packet (x = 1); (b) x = 0.6 (in the 
lower Alfven speed region); (c) x = 1.4 (in the higher Alfven speed region) for R = S = 10 4 , 
a = 2 and d = 0.1. 
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Figure 4: The shape of the packet at (a) t = 25s ; (b) t = 100s ; (c) t = 175s for R = S = 10 4 , 
a = 2 and d = 0.1. The dash-dotted curve is the initial gaussian packet. 
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Figure 5: The time dependencies of overall (a) kinetic energy; (b) magnetic energy and (c) the 
decay rate of total energy of the wave packet for R = S = 10 4 , a = 2 and d = 0.1. 
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(c) 
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Figure 7: Same as Fig. 4 for R = S = 10 s at (a) t = 500s ; (b) t = 1500s ; (c) t = 3000s. 
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Table 1: Numerical and analytical (based on HP83) values of A and B for some Reynolds and 
Lundquist numbers in the cases of S = R and S 3> R. Note that Analytic = 3 in HP83. 



R 


s 


A 

^numeric 




^analytic 




^numeric 


10 4 


10 4 


2.53 x 10" 


-3 


2.63 x 10- 


-7 


1.277 


10 5 


10 5 


1.73 x 10" 


-5 


2.63 x 10- 


-8 


2.001 


10 6 


10 6 


2.69 x 10" 


-8 


2.63 x 10- 


-9 


2.692 


10 7 


10 7 


1.89 x 10- 


-9 


2.63 x 10" 


-10 


2.914 


10 8 


10 8 


1.23 x 10- 


-10 


2.63 x 10" 


-11 


2.980 


10 9 


10 9 


1.08 x 10" 


-11 


2.63 x 10- 


-12 


2.992 


10 10 


10 10 


1.05 x 10- 


-12 


2.63 x 10" 


-13 


2.995 


10 7 


10 12 


7.70 x 10- 


-10 


1.32 x 10" 


-10 


2.946 


10 8 


10 12 


5.84 x 10- 


-11 


1.32 x 10" 


-11 


2.987 


10 9 


10 12 


5.36 x 10- 


-12 


1.32 x 10- 


-12 


2.993 
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